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Abstract. The stability properties of simple element choices for the mixed 
formulation of the Laplacian are investigated numerically. The clement choices 
studied use vector Lagrange elements, i.e., the space of continuous piecewise 
polynomial vector fields of degree at most r, for the vector variable, and the 
divergence of this space, which consists of discontinuous piecewise polynomials 
of one degree lower, for the scalar variable. For polynomial degrees r equal 2 
or 3, this pair of spaces was found to be stable for all mesh families tested. In 
particular, it is stable on diagonal mesh families, in contrast to its behaviour 
for the Stokes equations. For degree r equal 1, stability holds for some meshes, 
but not for others. Additionally, convergence was observed precisely for the 
methods that were observed to be stable. However, it seems that optimal order 
L 2 estimates for the vector variable, known to hold for r > 3, do not hold for 
lower degrees. 



1. Introduction 

In this note, we consider approximations of the mixed Laplace equations with 
Dirichlet boundary conditions: Given a source g, find the velocity u and the pressure 
p such that 

u — gradp = 0, div u — g in fi, p = on dQ, 
for a domain f2 C M 2 with boundary <9f2. The equations offer the classical weak 
formulation: Find a square integrable vector field with square integrable divergence 
u G if (div, Q) and a square integrable function p G L 2 (il) such that 



(1.1) u ■ v + q div u + p div v — I g q 

Jn Ja 

for all v G if (div,f2) and q G L 2 (il). The above formulation can be discrctized 

using a pair of finite dimensional spaces Vh C if (div, f2), Qh C L 2 (il), yielding 

discrete approximations uy L G Vh and Ph & Qh satisfying ( |1.1[ ) for all v G Vh and 

q G Q h - 

As is well-known, the spaces Vh and Qh must satisfy certain stability, or com- 
patibility, conditions for the discretization to be well-behaved 7 . More precisely, 
there must exist positive constants a and (3 such that for any h, 

(li, f) 

(1.2a) < a < ah — inf sup j. — r. — '-rr—r, — , 

uEZ h veZh |M!div|M|div 

(1.2b) 0<(3</3t= inf sup (di ™ ,<z) 



q£Qh veVh IMIdivlkllo' 
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Here, || • ||div and || • ||o denote the norms on H(div, SI) and L 2 (S7), (•, •) is the 
L 2 (S1) inner product and 

(1-3) Z h = {veV h \(divv,q)=0 VqeQ h }. 

The two conditions will be referred to as the Brezzi coercivity and the Brezzi inf-sup 



condition for the mixed Laplacian. The classical conforming discretizations of ( 1.1 ) 
rely on the finite element families of Raviart and Thomas [18 a or Brezzi, Douglas 
and Marini [8 j for the space Vh C iJ(div) in order to satisfy these conditions. 

In this note, we shall consider the Lagrange vector element spaces, that is, con- 
tinuous piecewise polynomial vector fields defined relative to a triangulation 7/j, 
for the space Vh- This is motivated by the following reasons. First, these spaces 
are fairly inexpensive, simple to implement and post-process and in frequent use 
for other purposes. Second, such pairs would allow continuous approximations of 
the velocity variable, or when viewed in connection with linear elasticity, lay the 
ground for continuous approximations of the stress tensor. Moreover, in the recent 
years, there has been an interest in mixed finite element discretizations that are 



both stable for (1.1 1 and for the Stokes equations: 



(1.4) / grad u : grad v + q div u + p div v = / fv 

Jn Jn 

for all v E i/ 1 (f2;V) such that J Q v = and all q E L 2 (S7). The search for con- 
forming such discretizations is complicated by the fact that the existing, stable 
discretizations of ( |1.1[ ) are such that Vh (£ i? 1 (fi;V). On the other hand, the 
existing stable discretizations of (1.4 1 are typically unstable for ( |1.1[ ) [13] . The ex- 

lcfnn/in r~vF cfol-iln r\ i cvrnf wall An c \/, v f) , i"v+ /1 1 I onciVi fVio^ T/, r~ TJ*-{\}'\I\ Kor>Ar 



istence of stable discretizations Vh x Qh of (1.1) such that Vh C i? 1 (SI; V) becomes 
a natural separate question. Unfortunately, there are no known such finite clement 
discretizations that are stable for any admissible triangulation family {7h}- In this 
note, we aim to numerically examine cases where a reduced stability property may 
be identified. In this sense, the investigations here are in the spirit of the work of 
Chapelle and Bathe [TU] and Qin [IT] . 

For a family of conforming discretizations {Vh x Qh}h of ( |1.1| such that div Vh C 



Qh for each h, the condition (1.2a) is trivial. The stability conditions thus reduce to 
the condition (1.2b I, namely the question of bounded Brezzi inf-sup constant /3^ IV - 



On the other hand, recall that for the Stokes formulation (1.4 1, the corresponding 
Brezzi coercivity condition is trivial by the Poincare inequality. Hence, for any 
family of conforming discretizations, the stability conditions for Stokes reduce to 
that of a uniform bound for the Brezzi inf-sup constant (3^. Here, 

(1.5) fa = inf sup . 

q£Qh veVh IMIilMlo 

when V h C iT^tyV), Q h C L 2 (S1) and || • ||i denotes the norm on _ff 1 (f2). Further, 



such a bound immediately gives (1.2b I since (3\ < /3jJ lv by definition. Hence, if 
divVh C Qh, stability for Stokes immediately gives stability for the mixed Lapla- 
cian. 

The conditions divV/j C Q h and Vh C i? 1 (Sl;V) are clearly satisfied by the 
element pairs consisting of continuous piecewise polynomial vector fields of degree 
less than or equal to r and discontinuous piecewise polynomials of degree r — 1, 
for r = 1,2, ... . This family could be viewed as an attractive family of elements 
for both the Stokes equations and the mixed Laplacian. However, the Brezzi inf- 
sup constant(s) will not be bounded for all r. For r > 4, Scott and Vogelius 
demonstrated that these finite element spaces will be stable for the Stokes equations 
on triangulations that have no nearly singular vertices, that is, triangulations that 
are not singular in the appropriate sense [19] , The lower order cases, 1 < r < 
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3, were studied carefully by Qin, concluding that the elements are not stable in 
general [17]. However, they are stable for some specific families of triangulations, 
and can be stabilized by removal of spurious pressure modes on some other classes 



of triangulations. (The space Mh of spurious pressure modes is defined in (2.10) 
below.) In general, the stability of finite element spaces for the Stokes equations 
has been extensively investigated. In addition to the previous references, surveys 
are presented in [3j |9]. However, to our knowledge, a careful study of the lower 
order cases has not been conducted for the mixed Laplacian. 

As the stability for the mixed Laplacian is a weaker requirement when div Vh Q 
Qh, there may be a greater class of triangulations for which the elements form a 
stable discretization. In fact, this is known to be true. One example is provided 
by the pairing of continuous piecewise linear vector fields Vh and the subspace of 
discontinuous piecewise constants Qh such that Qh — div Vh on crisscross triangu- 
lations of the unit square. Qin proved that there does not exist a (3 > such that 
Pi > (3 for any h [TTl Lemma 7.3.2]. On the other hand, Boffi et al. proved that 
such a bound does exists for /3^ lv [3] . We shall present numerical evidence suggest- 
ing that there is a range of triangulations for which this phenomenon occurs. The 
main results are summarized below. 

Spurious modes: For r = 2,3 and for all triangulations tested, the dimen- 
sion of the space Mh of spurious modes is equal to the number of interior 
singular vertices a. However, for r = 1 and one of the triangulation families 
studied ("Flipped", which is defined in Figure [l}, dimjV/j is strictly greater 
than a . 

Stability: For all triangulations we have tested, the method seems at least 
reduced stable (i.e., stable after removal of spurious modes, if any), for 
r — 2,3. This is in contrast to the situation for the Stokes equations, where 
for some triangulations, such as the diagonal triangulation, the method is 
not reduced stable for r = 2,3, while for other triangulations, it is. For 
r = 1, reduced stability holds for some triangulations, but fails for others, 
including the diagonal triangulation. 
Convergence: We also studied convergence of the method on diagonal tri- 
angulations. For such meshes, the method was observed to be stable for 
r > 1, but unstable for r = 1. Theory predicts optimal convergence of p 
in L 2 (il) and u in iJ(div,fi) for a stable method and this is in fact what 
was observed. Such optimal convergence holds for r = 2, 3, 4, but not in 
the apparently unstable case r = 1. In the case r > 4, it is known that u 
converges at one order higher in L 2 (il) than in iJ(div, fi). No such increase 
of order was observed for r < 4. 
The note is organized as follows. We introduce further notation and summarize 
some key points of the theory of mixed finite element methods in Section[2j Further, 
we derive some eigenvalue problems associated with the stability conditions and 
give a characterization of the Brezzi inf-sup constant for the mixed Laplacian /3jJ lv 
in Section [3j These eigenvalue problems applied to the Stokes equations were 
also stated by Malkus [13] and (in part) by Qin [T7] and provide a foundation for 
numerical investigations of the Brezzi stability conditions. Based on these general 
results, Section [4] is devoted to the study of continuous piecewise polynomials in 
two dimensions for the velocity and discontinuous piecewise polynomials for the 
pressure. 

2. Notation and preliminaries 



The notion of reduced stability of families of mixed finite element spaces is a key 
point in this note. In order to make this notion precise, this preliminary section 
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aims to introduce notation and summarize the stability notions for finite clement 
discretizations of abstract saddle point problems. 

If V is an inner product space, we denote the dual space by V* , the inner product 
on V by (• , • )y and the induced norm by || • ||y. Let f2 be an open and bounded 
domain in W 1 with boundary dfl. We let H m {Q), for m = 0, 1, ... , denote the 
standard Sobolev spaces of square integrable functions with m weak derivatives 
and denote their norm by || • \\ m . Accordingly if°(f2) = L 2 (il). The space of 
polynomials of degree r on ft is denoted 7-V(f2). The space of vectors in R d is 
denoted V and in general, X(Q; V) denotes the space of vector fields on for which 
each component is in X(fi). For brevity however, the space of vector fields in 
L 2 (fl; V) with square integrable divergence is written _ff(div,J7) with norm |j • ||<jiv 
and semi-norm | • |div = II div • 1 1 o - The subscripts and the reference to the domain 

will be omitted when considered superfluous. 

Let 7/,, denote an admissible simplicial tessellation of f2, h measuring the mesh 
size of the tessellation. We shall frequently refer to spaces of piecewise polynomials 
defined relative to such, and label the spaces of continuous, and discontinuous, 
piecewise polynomials of degree less than or equal to r as follows. 

V c r =V c r {T h ) = {p£H 1 {tt)\p\ K eVr(K) VKeT h } r=l,2,..., 
V r = V r (T h ) = {peL 2 (n)\p\ K eV r (K) VKeT h } r = 0,l,.... 

The classical abstract saddle point problem reads as follows [71 for given 
Hubert spaces V and Q and data (/, g) £ V* x Q* , find (u,p) £ V x Q satisfying 

(2.1) a(u,v) + b(v,p) + b{u,q)={f,v) + {g,q) V(v,q)£VxQ, 

where a and b are assumed to be continuous, bilinear forms on V x V and V x Q, 
respectively. We shall assume here and throughout that a is symmetric. Follow- 



ing [2], there exists a unique solution (u,p) of (2.1 1, if and only if the continuous 
Babuska inf-sup constant 

, > a(u,v) + b(v,p)+b(u,q) 

(2.2) 7= mf sup — 

0?(v,q) \\\ u ,P)\\vxQ\\(V,q)\\vxQ 

is positive. By [7], this holds if and only if the continuous Brezzi coercivity and 
Brezzi inf-sup constants are positive. These are defined as 

/ r\ q \ • t a(u,v) 

(2.3) a = mf sup 



(2.4) P= inf sup 



o^uez -£ veZ ||u||y||u||y ' 
b(v,q) 



o#9eQo^ev IMIv||<31|q 



respectively, where Z — {v £ V | b(v, q) = Vg £ Q}. 

Given finite dimensional spaces Vh C V and Qh C Q, defined relative to a 
tessellation 7^ of f2, the Galerkin discretization of (2.1l takes the form: Find 
(u h ,p h ) £ Vh x Q h satisfying 



(2.5) a(u h ,v) + b(v,p h ) + b(u h ,q) = (f,v) + (g,q) V (v, q) £ V h X Q h . 



STABILITY OF LAGRANGE ELEMENTS FOR THE MIXED LAPLACIAN 



5 



On the discrete level, the Babuska inf-sup, Brezzi coercivity and Brezzi inf-sup 
constants are defined as 



(2.6) 7/j := inf sup 



a(u, v) + b(v,p) + b(u, q) 



0^(u, P )eV h xQ h ^ {Vyq)eVhXQh \\(u,p)\\ VxQ \\(v,q)\\ v> 

fo-7\ ■ t a(u,v) 

(2.7) a h := mf sup , 

o^uez h0 ^ veZh \\u\\v\\v\\v 

(2.8) (3h ■= inf sup f — , 

o^q£Q h0 ^ veVh \\v\\vH\q 

where 

(2.9) Z h = {veV h \b(v,q) = VqeQ h }. 



For given Vh x Qh, there exists a unique solution of (2.5 1 if and only if and fih 
(or equivalently 7^) are positive. Furthermore, for a family of discretization spaces 
Vh x Qh parameterized by h, if ah and f3h are uniformly bounded from below, then 
one obtains the quasi-optimal approximation estimate [7]: 



\\u-u h \\v + \\p-PhWq < C inf \\u - v\\ v + inf \\p-q\\ Q 

\veV h qeQ h 

with C depending only on the bounds for ah and f3h and bounds on the bilinear 
forms a and b. The uniform boundcdncss condition motivates the notion of stability 
for pairs of finite clement spaces. 

Definition 2.1 (Stable discretization). A family of finite element discretizations 
{Vh x Qh}h is stable in V x Q if the Brezzi coercivity and inf-sup constants ah and 
(3h (or equivalently the Babuska inf-sup constant jh) o,re bounded from below by a 
positive constant independent of h. 

In accordance with standard terminology, we say that {Vh x Qh}h satisfies the Brezzi 
coercivity or inf-sup conditions if ah or fth, respectively, are uniformly bounded from 
below. 

There arc families of discretizations that are not stable in the sense defined above, 
but have a reduced stability property. More precisely, for a pair Vh x Qh consider 
the space of spurious modes Mu Q Qh'- 

(2.10) M h = {q£Q h \b(v,q)=0 Vt> S Vh}. 

For a stable discretization, A//, contains only the zero element. Indeed, (3h = if 
and only if Mh contains non-zero elements. If Mh is non-trivial, it is natural to 
consider the reduced space Af^, the orthogonal complement of Afh in Qh, in place 
of Qh- This motivates the definition of the reduced Brezzi inf-sup constant, relating 
to the stability of Vh x Af^: 

(2.11) h = inf sup h ^°' q ' 



and the following definition of reduced stable. By definition, f3h ^ 0. 

Definition 2.2 (Reduced stable discretization). A family of discretizations {Vh x 
Qh}h is reduced stable in V x Q if the Brezzi coercivity constant ah and the reduced 



Brezzi inf-sup constant (3h, defined by (2.7| and (2.111, are bounded from below by 
a positive constant independent of h. 
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3. Eigenvalue problems related to the Babuska-Brezzi constants 



For a given set of discrete spaces, the Babuska and Brezzi constants denned 
by ( 2.6 1 — ( 2.8) can be computed by means of eigenvalue problems. The form and 



properties of the eigenvalue problem associated with the Brezzi inf-sup constant 
for the Stokes equations were discussed by Qin in [T7]. Since also the Brezzi co- 
ercivity constant plays a role for the mixed Laplacian, we begin this section by 
deriving how the Brezzi coercivity constant can be computed by similar eigenvalue 
problems. Actually, in our application in Section [4] the Brezzi coercivity condition 
will be automatic, however, we discuss it here in the abstract case, for the sake of 
completeness. These eigenvalue problems were also stated, and carefully analysed 
from an algebraic view-point, by Malkus |13j in connection with the displacement- 
pressure formulation of the linear elasticity equations. We continue by observing 
that the continuous Brezzi inf-sup constant can be naturally associated with the 
smallest eigenvalue of the Laplacian itself. 

3.1. Eigenvalue problems for the discrete Babuska-Brezzi constants. Let 

Vh C V and Qh C Q be given finite dimensional spaces as before. It follows easily 
from the definition that the Babuska inf-sup constant jh — |A m in| when A m ; n is the 
smallest (in modulus) eigenvalue of the following generalized eigenvalue problem: 
Find A G K, ^ (u,p) G Vh x Q h satisfying 

(3.1) a(u, v) + b(v,p) + b(u, q) = X((u,v) v + (p,q)q) V (v,q) G Vh X Q h . 

The following lemma identifies an eigenvalue problem associated with the Brezzi 
inf-sup constant. 

Lemma 3.1 (Qin [IT , Lemma 5.1.1 - 5.1.2]). Let A m j n be the smallest eigenvalue 
of the following generalized eigenvalue problem: Find A £ M, 0^ (u,p) € Vh X Qh 
satisfying 

(3.2) {u,v) v + b(v,p)+b(u,q) = -X(p,q) Q , V («, q) G V h x Q h . 
Then, A > and for (3h defined by (2.8 1, j3h = vAmin- 



It can also be shown that the reduced Brezzi inf-sup constant (3h equals the square- 



root of the smallest non-zero eigenvalue of (3.2 1 [17, Theorem 5.1.1]. 

The Babuska and Brezzi inf-sup constants are thus easily computed, given bases 



for the spaces Vh and Qh- As for (3.1 1, it is easily seen that the Brezzi coercivity 
constant a/, = |A m j n | where A m i n is the smallest (in modulus) eigenvalue of the 
eigenvalue problem: Find A £ i and ^ u G Zh such that 

(3.3) a(u, v) = X(u, v)v V v G Zh- 

However, a basis for Zh is usually not readily available, thus hindering the actual 



computation of the eigenvalues of (3.3 1. Instead, the above eigenvalue problem over 



Zh can be extended to a generalized eigenvalue problem over Vh x Qh'- Find A G 
and 7^ (u,p) G Vh x Qh such that 



(3.4) 



a(u, v) + b(v,p) + b(u, q) = X(u, v) v V (v, q) G V h X Q h . 



The following lemma establishes the equivalence between (3.3 1 and (3.4) 



Lemma 3.2. If (A, it) is an eigenpair of (3.3), there exists a p G Qh such that 
(A, (u,p)) is an eigenpair of ( 3.4 1 . Conversely, if {X, {u,p)) is an eigenpair of (3.4) 
and «/0, then u G Zh and (A, u) is an eigenpair of (3.3 1. For 0^p£ Afh and 
any scalar X, (A, (0,p)) is an eigenpair of (3.4), and these are the only eigenpairs 
of (pL4) with u = 0. 
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-> Bh{Vh) is an isomorphism, 
\/q&B h {V h ). 



Proof. Let (A, u) be an eigenpair of (3.3 1. Define Bh : Vh — * such that 
(B h v,q) Q = b(v,q) for all g G Q^. Since : 
p G B h (Vh) C <5^ is well-defined by 

(Pj g>Q = A(u, B~ 1 q) v - a(u, B~ 1 q) 

Then, for any v G Z^, p satisfies 

H V ,P) = ( B hV,p)Q = X(u,v) v - a(u,v). 

Further, by definition b(v,p) — for any v G 2%. Hence, by the assumption that 
(A, w) is an eigenpair of ( |3.3| ), (A, (u,p)) satisfies (3.4|. The converse statement is 
obvious. Finally, letting u — in ( |3.4[ ), we see that (A, (0,p) satisfies (3.4 1 if and 
only if p G Afh, but for any Ael. □ 

Note that, as a consequence of the last observation in Lemma [3.2| if Afh is non- 



trivial, the generalized eigenvalue problem (3.4 1 is computationally not well-posed 
since any scalar A is an eigenvalue. 

In the subsequent section, we shall numerically investigate the stability of families 
of finite element discretizations Vh X Qh such that div Vh C Qh for the mixed 
Laplacian, using the eigenvalue problem (3.2) in terms of standard bases for the 
spaces Vh and Qh- The eigenvalue problem (3.4 1 does not enter, but would, were 
we to investigate discretizations where div Vh % Qh- 

3.2. A characterization of the mixed Laplacian Brezzi inf-sup constant. 

We now turn from the general setting to consider the H(div) x L 2 formulation of 
the mixed Laplacian (1.1). In Lemma 3.3 below, we show that the Brezzi inf-sup 
constant can be identified with the smallest eigenvalue of the negative Laplacian. 
Consequently, if a discretization family {Vh X Qh}h guarantees eigenvalue conver- 
gence for the mixed Laplace eigenvalue problem and is such that div Vh Q Qh , the 
Brezzi inf-sup constant of the discretization will converge to the continuous Brezzi 
inf-sup constant. 

Lemma 3.3. Let V C iJ(div, f2) and Q C L 2 (fl) be such that that divV C Q. 
Consider the Brezzi inf-sup eigenvalue problem (3.2) applied to (JlTJ): 

(3.5) (u, v)div + (divv,p) + (divu,g) = -A(p, g) y(v,q)eVxQ. 
Consider also the mixed Laplace eigenvalue problem: 

(3.6) (u,v) + (divv,p) + <divn,g) = -X{p,q) V(v,q) eV x Q. 



Then, (X,(u,p)) is an eigenpair of (3.51 if and only if (A, (u,p)) is an eigenpair 
of (3.6 1 where A = A(l — A) -1 , u = u and p = (1 — X)p. Moreover, in this case, 
< A < 1 and A > 0. Also, p^O.p^O. 



Proof. Assume that (X,(u,p)) is an eigenpair of (3.51. First, note that A ^ 1. 
Letting A = 1, v — u and q = — div u in (3.51, implies that u — 0. Further, 
v = and q = p gives that p = 0. Hence, A = 1 is only associated with the zero 
solution, which by definition, cannot form an eigenpair. Next, note that p ^ 0, 
since otherwise implies that u — 0, which again is impossible. By the assumption 
divV C Q, 

(3.7) (div u, div w) = — X(p, divv) Vu G V. 

Taking v = u in (3.7l and letting v — u and g = (A — l)p in (3.51 show that 
||u|| 2 = A(l - A)||p|p~So, < A < 1. The combination of {3J| and (p3J5| gives 

(u,v) + (divw,g) + (1 - X)(dxvv,p) = -X(p,q). 



Finally, letting u = u, p = (1 — X)p and A = A(l — A) 1 , 
solves (3.6). The converse holds by similar arguments. 



jives that (A, (u,p)) 
□ 
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The equivalence demonstrated in the lemma above affords a simple character- 
ization of the Brezzi inf-sup constant for the mixed Laplacian. The eigenvalue 
problem (3.5 1 with V = -ff(div,f2) and Q = L 2 (f2) is the eigenvalue prob lem as- 
sociated with the continuous Brezzi inf-sup constant /3 dlv , cf. Lemma |3.l| Hence, 
/3 dlv is the square-root of the smallest eigenvalue of (|3.5| . On the other hand, the 



eigenvalue problem (3.6 1 is a mixed weak formulation of the standard eigenvalue 



problem for the negative Laplacian with Dirichlet boundary conditions, given in 
strong form below: 

(3.8) - Ap = Xpinn, p = 0ondCl. 

Thus, if p is the smallest eigenvalue of (3.8 1, /3 dlv = \J /i(l + 



Remark. An alternative eigenvalue problem arises from noting that, under the as- 
sumption div V C Q, (3.6) with v = implies that Xp = — divu. Hence, if 



(X,(u,p)) solves (3.6 ), then cither A = 0, u = and p _L div V, or A > 0, u ^ 0, and 
(A, u) will be an eigenpair of the problem: 

(div u, div v) = X(u, v) V v € V. 
This eigenproblem was studied in [5]. 



Now, consider a stable discretization family Vh x Qh of (1.1 1 such that div Vh Q 
Qh, with Brezzi inf-sup constants /3^ lv - Let p^ denote the smallest eigenvalue ap- 
proximation of (3.6 1 by Vh X Qh- As a consequence of the preceding considerations, 
if fj,h — ► /i, then /3JJ 1V — > /3 dlv . In other words, if the discretization family is stable, 
satisfies div 14 C Q h , and gives eigenvalue convergence, then the Brezzi inf-sup 
constant will converge to the continuous Brezzi inf-sup constant. Note however, 
that the discrete stability conditions are not sufficient to ensure the convergence of 



approximations to the eigenvalue problem (3.6 1 [UH]. 

Mixed finite element discretizations of (1.1 1 based on the Raviart-Thomas 



and Brezzi-Douglas-Marini [5] families of H (div) conforming elements are known 
to give eigenvalue convergence, and hence /3j] lv — > /3 dlv . Some cases where /3jJ lv 
seems to be uniformly bounded in h, but /? ; dlv f3 dlv are exemplified in the subse- 
quent section. Finally, note that if is the unit square: f2 = (0, l) 2 , the smallest 



eigenvalue of (3.8 1 is 2ir 2 and so 



(3.9) 



div 



2tt 2 



1 + 2tt 2 



0.975593. 



4. Lower order Lagrange elements for the mixed Laplacian 

From here on, we restrict our attention to finite clement discretizations of the 
mixed Laplacian (1.1) on a polygonal domain n C I 2 . The primary aim is to 
examine the stability, or reduced stability, and convergence properties of Lagrange 
elements, that is, continuous piecewise polynomials for the vector variable and 
discontinuous piecewise polynomials for the scalar variable: 



(4.1) 

for r = 1,2, 



V h xQ h = V c r {T h ;N) x V r -i{T h ), 



Although the Brezzi conditions are in general not satisfied for 
these discretizations, stability or reduced stability may be identified on families 
of structured triangulations. The pair (4.1) is clearly such that divV/j C Q h . 



Therefore, the stability of the discretization relics on a uniform bound for the 
Brezzi inf-sup constant only. Further, a uniform lower bound on the Brezzi inf-sup 
constant for the Stokes equations induces the corresponding bound for the mixed 
Laplacian. Hence, the results on the reduced stability of this clement pair for the 
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Stokes equations can be directly applied to the mixed Laplacian. In the following, 
new numerical evidence is presented and compared to the known results. 

The stability of the V^(V) x V r -i family of elements, for both the Stokes equa- 
tions and the mixed Laplacian, depends on the polynomial degree r and the struc- 
ture of the triangulation 7/j. For triangulations that have interior singular vertices, 



the space of spurious modes Afh, defined by (2.101 applied to (1.1 1, will be non- 
trivial. Here, an interior vertex is labeled singular if the edges meeting at that 
vertex fall on two straight lines. Let x denote an interior singular vertex and let oj x 
be the star of x. For any r > 1, there exists ape Afh such that p is supported in 
w x |16l 115] . Consequently, letting a denote the number of interior singular vertices 
of a triangulation, dimA/h > o. Scott and Vogelius showed, for r > 4 that if there 
are no interior singular ver tices, then dimA/k = and so [if v > (3\ > 0. More- 
over, they proved that for a family of meshes without interior singular vertices, (3\ 
remains bounded above zero as long as the meshes do not tend to singularity as 
h — > 0. For the precise statement and more details, see [T!5] or Section 12.6]. 

As we shall see below, for r < 4, the space of spurious modes may be non-trivial 
even when there are no singular vertices. Further, for the Stokes equations, more 
restrictive conditions than the above must be placed on the triangulations in order 
to obtain a uniform bound for the Stokes Brezzi inf-sup constant [17 . The stability 
properties of these lower order discretizations for the mixed Laplacian is the main 
question of interest in the following. 

Remark. We shall not consider the pairing P£(V) x V s except for s = r — 1. This 
choice is easily motivated. First, a dimension count shows that the pairing of 
continuous piecewise polynomial vector fields with discontinuous piecewise polyno- 
mials of the same or higher degree must have a non-trivial space of spurious modes. 
Second, although the Brezzi inf-sup constant is uniformly bounded for the pairs 
V£(y) x Vr-2, t = 2, 3, . . . , the Brezzi coercivity constant for the mixed Laplacian 
is not uniformly bounded, and thus stability fails. 

4.1. Stability. In the spirit of [I7J Section 5], we aim to numerically investigate 
the stability of V£(V) x V r -\ for r = 1, 2, 3 on certain families of structured triangu- 
lations of the unit square. The triangulation patterns considered are illustrated and 
labeled in Figure [l] For n even, annxn triangulation of each family is constructed 
by first partitioning the domain into n x n squares, and subsequently dividing each 
block of 2 x 2 squares into triangles by the respective patterns. For instance, annxn 
diagonal triangulation is formed by dividing the unit square into n x n subsquares, 
and dividing each subsquare into triangles by the positive diagonal. Throughout, 
we identify h = 1/n and assume that n > 2. Observe that the diagonal, flipped, and 
zigzag triangulations contain no interior singular vertices, while the crisscross and 
the Union Jack triangulation contain n 2 and n{n — 2) /2 interior singular vertices, 
respectively. This is summarized in the first row of Table [T] 

Recall that dimA/fc > er for r > 1 and equality holds for r > 4. Qin proved 
that equality holds for 1 < r < 3 in the case of the diagonal and the crisscross 
meshes and numerically observed equality for the flipped mesh for r = 2 [17] . Our 
own experiments show that equality holds for the zigzag and Union Jack meshes 
for 1 < r < 3. Equality also holds for the flipped mesh when r = 2,3, but not for 
r = 1. These results are summarized in the second row of Table [TJ 

We continue by studying the behaviour of the Brezzi inf-sup constants on the 
above triangulations. The cases r = 2, 3 are considered first, but we will return to 
the case r = 1 below. For the Stokes equations, it is known that the diagonal and 
crisscross triangulation families exhibit very different behaviour for r — 2,3 [TT] . 
Namely, although there are non-trivial spurious modes on the crisscross triangula- 
tion family, the reduced Brezzi inf-sup constant is uniformly bounded. In contrast, 
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(a) Diagonal (b) Flipped (c) Zigzag 



(d) Crisscross (c) Union Jack 



Figure 1. Structured 2x2 triangulations of the unit square. 





Diagonal 
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Zigzag 


Crisscross 


Union Jack 


a 













o 

n 


n{n- 2)/2 


dim Afh 





Un/2 


-l) 2 r=l 





n 2 


n(n-2)/2 






1° 


r = 2,3 







Table 1 . The number of interior singular vertices a and the di- 
mension of the space of spurious modes, dim Afh, for labeled fami- 
lies ofnxn triangulations of the unit square, cf. Figure [T] For the 
flipped, zigzag and Union Jack meshes, dimA/h, is conjectural. 



for the diagonal family, the Brezzi inf-sup constant decays as approximately 0(h). 
As the discretization is reduced stable for the Stokes equations on crisscross tri- 
angulations, it is also reduced stable for the mixed Laplacian. A natural question 
becomes whether the lack of stability on diagonal triangulations for the Stokes 
equations is also present for the mixed Laplacian. 

In view of Lemma |3.1[ we shall make an attempt at answering this question 
through a set of numerical experiments. For a given r and a given 7^, the small- 



est, and smallest non-zero, eigenvalue of (3.5 1 for V — V£(Th,V), Q = V r -i{Th) 
give the Brezzi inf-sup and reduced Brezzi inf-sup constant. These eigenvalues for 
the triangulation families considered, computed using LAPACK, SLEPc [TT] and 
DOLFIN [121 . are given for r = 1, 2, 3 in Tables |J [3] and g For the purpose of iden- 
tifying spurious modes, eigenvalues below a threshold of 1CP 4 have been tabulated 
to zercO 

For the diagonal meshes, the numerical experiments indicate that in contrast 
to Stokes, the mixed Laplacian Brezzi inf-sup constants are bounded from below 
for both r = 2,3. For the flipped and zigzag meshes, experiments give similar 

^Had a smaller threshold been chosen, some of the zero eigenvalues associated to interior 
singular vertices would have been missed for r = 2 on the Union Jack mesh of size n = 6. 
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(dimMO 



Union Jack 



4 


0.847171 


0.791967 


0.945496 (1) 


0.976985 


(4) 


6 


0.716677 


0.626865 


0.945619 (4) 


0.976271 


(12) 


8 


0.605576 


0.505968 


0.947850 (9) 


0.975985 


(24) 


10 


0.517707 


0.420180 


0.946138 (16) 


0.975847 


(40) 


12 


0.449060 


0.357720 


0.944833 (25) 


0.975770 


(60) 


14 


0.394963 


0.310731 


0.943880 (36) 


0.975724 


(84) 


16 


0.351684 


0.274303 


0.943142 (49) 


0.975693 (112) 



Table 2. The mixed Laplacian (reduced) Brczzi inf-sup constant 
for Vi(Th;V) x Vo{Tk) on labeled structured families of triangula- 
tions Th . The dimension of the space of spurious modes in paren- 
thesis if non-trivial. 
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(dimMO 
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Diagonal Zigzag Flipped Union Jack 


4 
6 
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10 
12 
14 


0.975627 0.955956 0.943790 
0.975600 0.952460 0.940480 
0.975595 0.951384 0.938717 
0.975594 0.950906 0.937684 
0.975594 0.950638 0.936992 
0.975593 0.950458 


0.975628 (4) 
0.975603 (12) 
0.975595 (24) 
0.975594 (40) 
0.975593 (60) 



Table 3. The mixed Laplacian (reduced) Brczzi inf-sup constant 
for V^T^V) x Vi(Th) on labeled structured families of triangula- 
tions Th. The dimension of the space of spurious modes in paren- 
thesis if non-trivial. 







(dimMO 


n 


Diagonal 


Zigzag 


Flipped 


Union Jack 


4 


0.972244 
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0.975594 


0.975594 (4) 


6 


0.967304 
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0.964845 


0.975593 


0.975593 


0.975593 (24) 
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0.963412 








12 


0.962484 









Table 4. The mixed Laplacian (reduced) Brczzi inf-sup constant 
for V^(Th\ V) x ViiTh) on labeled structured families of triangula- 
tions Th . The dimension of the space of spurious modes in paren- 
thesis if non-trivial. 
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results. Neither exhibits any spurious modes. Moreover, while the Brezzi inf-sup 
constant decays approximately as O(h) for the Stokes equations [17], it appears 
to be uniformly bounded for the mixed Laplacian. For the Union Jack family, the 
same experiment gives n(n — 2)/2 spurious modes, but the reduced Brezzi inf-sup 
constant again seems to be uniformly bounded. In summary for r = 2,3, the 
V£(V) x V r -i elements appear to be at least reduced stable for all the families 
considered. 

With the discussion in Section [3. 2 1 in mind, we also note that the Brezzi inf-sup 



constant converges to the exact value, given by (3.9), for some, but not all, of these 
meshes. For r — 2, the Brezzi inf-sup constant seems to converge to the exact value 
on the diagonal meshes, but not on the flipped or the zigzag meshes. The situation 
is the opposite for r — 3. There, the Brezzi inf-sup constant seems to converge to 
the exact value on the zigzag and flipped meshes, but not for the diagonal meshes. 

The situation is different and more diverse in the lowest-order case: r = 1. Bom 
et al. proved that "P° (V) x Pq is in fact reduced stable for the mixed Laplacian on 
crisscross meshes [I]. It is not reduced stable for Stokes [T7]. However, the element 
pair does not seem to be stable on diagonal meshes. The values in the first column of 
Table[2]indicate that the Brezzi inf-sup constant decays approximately as O(h). The 
same is the case for the zigzag meshes. For the Union Jack meshes, the situation 
is similar to the crisscross case. That is, the number of singular modes match 
the number of interior singular vertices and the reduced Brezzi inf-sup constant 
appears to be bounded from below. Finally, the flipped meshes display a surprising 
behaviour. There seem to be (n/2 — l) 2 spurious modes, even though there are 
no singular vertices. This is the only case where we have observed dim Afh > a. 
However, the reduced Brezzi inf-sup constant appears to be uniformly bounded. 

4.2. Convergence. In the previous subsections, we have investigated the stability 
of the V^(Y) x V r -i elements. Now, we proceed to examine the convergence prop- 
erties of these elements on the diagonal meshes. Conjecturing that V^(W) x V r -i 
is stable on this mesh family for r > 2, in accordance with the numerical evidence 
presented above, the standard theory gives the error estimate 

(4.2) \\u-u h \\ Aiv + lb-Pfc||o < C7i r (||u|| r+ i + ||p|| r ). 

For r > 4, the L 2 error estimate for the velocity can be improved [BJ Theorem 
12.4.9], thus yielding: 

(4.3) \\u-u h \\ <Ch r+1 \\u\\ r+1 . 



In order to verify (4.2) and to see whether (4.3) appears to be attained for 
r = 2,3, we consider a standard smooth exact solution to the Laplacian with pure 
Dirichlct boundary conditions: 

(4.4) p(x,y) = sin(27ra;) sin(27ry), u = gradp, g = divu. 

The errors of the 7*r(V) x V r -i approximations for 1 < r < 4 on diagonal meshes 
can be examined in Figure[2j To compute the errors, both the source function g and 
the exact solutions it, p have been represented by sixth order piecewise polynomial 
interpolants Tlhg and TlhU, HhP, whereupon the errors have been calculated exactly 
(up to numerical precision). 

For r = 1, we observed the discretization to be unstable on diagonal meshes. As 
expected in this case, neither the pressure nor the velocity approximation seems to 
converge in the L 2 norm. This indicates that the estimates for the approximation 
error, based on the standard estimates and the decaying Brezzi inf-sup constant, 
cannot be improved. On the other hand, for r > 1, we observed the method to be 
stable. For r = 2, 3, 4, the orders of convergence in the H (div) norm of the velocity 
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Figure 2. The errors of V^(V) x Vr-i approximations for r = 
1,2,3,4 on diagonal meshes versus mesh number n. The errors 
have been normalized, that is, multiplied by the inverse of the 
error at the smallest mesh n — 4. 
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and the L 2 norm of the pressure approximations are indeed optimal, as predicted 



by (4.21 



The situation seems different for the convergence of the velocity approximation 



in the L 2 norm. For r > 4, a convergence rate of order r + 1 is predicted by (4.3 1 



This is also observed for r = 4 in Figure 2(c) On the other hand, for r = 2 and 3 



the rate of convergence appears to be of order r, and thus one order suboptimal, 



indicating that the estimate (4.3 1 does not hold for r < 4. 
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